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ABSTRACT 

We investigate several approaches for constructing Monte Carlo realizations of the 
merging history of virialized dark matter halos ("merger trees") using the extended 
Press-Schechter formalism. We describe several unsuccessful methods in order to illus- 
trate some of the difficult aspects of this problem. We develop a practical method that 
leads to the reconstruction of the mean quantities that can be derived from the Press- 
Schechter model. This method is convenient, computationally efficient, and works for 
any power spectrum or background cosmology. In addition, we investigate statistics 
that describe the distribution of the number of progenitors and their masses as a 
function of rcdshift. 
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1 INTRODUCTION 

In the standard picture of modern structure formation, 
small-amplitude Gaussian density fluctuations, which per- 
haps arose from quantum fluctuations and were amplified 
by a period of rapid inflation, become more overdense with 
respect to their surroundings as the universe expands. Even- 
tually the self-gravity acting on these regions becomes larger 



than ti e pressure of the expansion, and they collapse to 
form bound, virialized structures. In hierarchical models, 
such as the cold dark matter (CDM) family of models, the 
amplitude of the fluctuations decreases with increasing scale. 
Thus small mass objects form first, and are then incorpo- 
rated into larger structures as time progresses. These dense, 
gravitationally bound structures provide the environments 
where galaxies can form. Hierarchical structure formation 
thus gives a natural explanation for the very complex ob- 
served large scale structure of the Universe, i.e. clusters, 



universe with Gaussian random phase initial perturbations. 
The focus of the original Press-Schechter model was the 
derivation of the multiplicity function of non-linear objects 
("halos") as a function of redshift, i.e. the "mass function" 
or number density of halos of a given mass at a redshift z. 
This prediction has been tested quite extensively and found 
t o be in relatively good agreement with N-body simulations 
( j Efstathiou et al. 198S|; |Gelb fc Berts chinger 1994 |Lacey & 



supercl isters, filaments, etc. 



One way to study this process is with N-body simula- 
tions. However, numerical simulations have familiar draw- 
backs. They are computationally expensive, so it is difficult 
or impossible to explore a wide range of models or differ- 
ent realizations of the same model. In addition memory and 
time limitations make it impossible to attain the mass and 
force resolution required to simultaneously study objects 
from dwarf galaxies (~ 10 9 Af Q ) to clusters (~ 10 15 M Q ). 
Semi-analytic methods are therefore an important alterna- 
tive. 



The model developed by Press & Schechter (1974) pro- 
vides a simple but relatively effective framework for the de- 
scription of the mass history of particles in a hierarchical 



Cole 1994 



Ma 1996 



Gross et al. 1998). The Press-Schechter 



theory was extended to give the conditional probability that 
a particle in a halo of mass Mq at zq was in a halo of mass 
Mi at an earlier redshift z\ . leading to an expression for th e 



conditional mass function (Bower 1991 



Bond et al. 



1991|) 



The extended Press-Schechter formalism can also be manip- 
ulated to obtain expressions f or halo survival tim es, forma- 
tion times, and merger rates (Lacey fc Cole 1993, hereafter 



LC93), which have also been shown to agree r easonably well 
with the results from N-body simulations (Lacey & Cole 



1994| ). 

The computation of these mean quantities within the 
Press-Schechter model is straightforward. However, for cer- 
tain purposes one would like to go beyond this. In particular, 
the semi-analytic approach to modeling galaxy formation 



(cf. Kauffmann, White, & Guiderdoni 1993, Cole et al. 1994) 



attempts to describe the formation history of galaxies and 
gas within dark matter halos, including simplified hydrody- 
namics, star formation, supernova feedback, galaxy-galaxy 
merging, and stellar population synthesis. These models rely 
on the construction of a "merger tree" , which involves pre- 
dicting the masses of progenitor halos and the redshifts at 
which they merge to form larger halos. Galaxies initially 
form in their own halo and are traced as they are incorpo- 
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rated into larger halos, and eventually perhaps merge with 
other galaxies. A halo of a given mass may have a variety of 
merging histories, and the properties of galaxies that form 
within this halo presumably depend to some extent on the 
details of this history. 

Most of the previous work using semi-analytic models 
has focussed on reproducing or predicting mean quantities 
and qualitative trends. However, as observational data con- 
tinues to improve, one would like to be able to investigate 
whether the broader properties of the predicted distribu- 
tion of model galaxies are consistent with the observations. 
For example, there has already been some investigation of 
whether the scatter in t he observed Tully-Fisher relation 
(Eiscnstein fc Loeb 1996), and in the color-magnitu de and 



line-str engt h velocity-dispersion (Mg-a) relations ( Kauff- 
mann ^996J) in merger models is consistent with observa- 
tions. Before we can trust the models for evaluating these 
kinds of questions, we must ensure that the merger trees not 
only satisfy the mean properties readily predicted by the ex- 
tended Press-Schechter model, but ideally also the full dis- 
tribution function. As we shall see, accomplishing this goal 
is far from straightforward and has not been thoroughly in- 
vestigated in previous work on this subject. 

In this paper, we discuss some of the practical and the- 
oretical difficulties of using the extended Press-Schechter 
model to create Monte Carlo merger histories of dark matter 
halos. We mention some limitations of the previously pro- 
posed methods for the construction of merger trees, and the 
motivation for developing a new approach. We discuss sev- 
eral unsuccessful approaches, in an attempt to clarify some 
aspects of this problem as well as to prevent others from 
following the same dead ends. In addition, the formalism 
we present may be useful in future work on this subject, as 
the final approach that we present seems to be effective and 
convenient, but it is still not rigorous. A primary motiva- 
tion for embarking on this project was to develop a method 
that reproduces the full joint probability distribution, rather 
than just the mean. Our method is not guaranteed to do so, 
and further investigation of this question requires a com- 
parison with N-bod y simulations. This wil l be presented in 
a companion paper (Somerville et al. 1998). 

In we give a brief introduction to the Press-Schechter 
formalism. In ^ we summarize some of the previous meth- 
ods for creating merger trees. In we describe the source of 
some difficulties one encounters in attempting to use the ex- 
tended Press-Schechter formalism to construct merger trees. 
We develop a simple model of the joint probability distribu- 
tion in §|H|, and attempt to use it to construct merger trees. In 
§[| we present a completely different approach which even- 
tually leads to our successful method, described in § |6 . 3| . In 
§0, we investigate the distribution of progenitor number and 
mass given by our successful method. We summarize and 
conclude in §H. 

Readers who are only interested in the successful 
method may skip to 



2 THE PRESS-SCHECHTER FORMALISM 



The Press-Schechter model ( Press fc Schechter 1974 ) 



based on a combination of linear growth theory, spheri- 
cal collapse theory, and the properties of Gaussian random 



fields. Suppose that we have smoothed the initial density 
distribution on a scale R using some spherically symmetric 
window function Wm(t), where M(R) is the average mass 
contained within the window function. There are various 
possible choices fo r the form of the window function (cf. 
Lacey & Cole 1993), and the relation between M and R will 



clearly depend on this choice. We use a real-space top hat 
window function, W M (r) = Q(R- r)(47i\R 3 /3)"\ where 6 
is the Heaviside step function. In this case M = 4-k Po R 3 /3, 
where po is the mean mass density of the universe. The mass 
variance S(M) = a 2 (M) may be calculated from 



a(M) = — / P{k)W\kR)k 2 dk, 



(1) 



where P(k) is the mass power spectrum, and W(kR) is the 
Fourier transform of the real space top-hat: 



W(kR) = 



3[sin(fc.R) - kRcos(kR)] 

(Why ' 



(2) 



The "excursion set" derivation due to Bond et al. (1991) 
leads naturally to the extended Press-Schechter formalism 
that we will use extensively in this paper. The smoothed 
field 8(M) is a Gaussian random variable with zero mean 
and variance 5*. The value of 5 executes a random walk as the 
smoothing scale is changed. Adopting an ansatz similar to 
that of the original Press-Schechter model, we associate the 
fraction of matter in collapsed objects in the mass interval 
M, M + dM at time t with the fraction of trajectories that 
make their first upcrossing through the threshold ui = S c (t) 
in the interval S,S + dS. This may be translated to a mass 
interval through equation (|l|) . The halo multiplicity function 
(here in the notation of LC93) is then: 



f(S,u)dS = 



1 



2^/2 



exp 



~2S 



dS. 



(3) 



The conditional mass function, the fraction of the tra- 
jectories in halos with mass Mi at z\ that are in halos with 
mass Mo at zq (Mi < Mo, Zo < Zi) is 



1 (wi 



2^ (Si - So) 3 ? : 



■ exp 



(ui — cjq) 2 
'2(5i -5* ) 



dSi. (4) 



The probability that a halo of mass Mo at redshift zo had 
a progenitor in the mass range (Mi, Mi + dMi) is given by 
(LC93): 



dP 

dMi 



(Ml, 31 



M ,z )dMi 

Mo 
Mi 



/(Si, wi I So, wo) 



dS 



dM 



dMi , (5) 



where the factor Mo /Mi converts the counting from mass 
weighting to number weighting. 

All of the results presented in this paper have been cal- 
culated for an Q — 1 universe with Ho = 50 km s _1 Mpc~ . 
The power spectrum is obtained from the fitting formula of 
Bardeen et al. (19861) with T = 0.21 and o s = 0-6. This is 



the tCDM model of Efstathiou, Bond, fc White (1992), and 



has been chosen because the slope of the power spectrum on 
galaxy scales is consistent with observations. However, our 
results are equally valid for any assumed power spectrum or 
cosmology. The one exception that we know of is the case of 
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a universe with a "hot" dark matter component, such as a 
massive neutrino (CHDM or MDM type models). The stan- 
dard extended Press-Schechter formalism does not properly 
treat the evolution due to the changing free-streaming length 
of the neutrino in such models. All of the expressions given 
are valid for a general cosmology unless otherwise noted. 



modeling, in which one would like to follow individual galax- 
ies with fairly fine time resolution. LC93 propose a general- 
ization of the block model which removes this condition, but 



we show in Section 3.1 that this method produces halo mass 



distributions that are severely discrepant with the Press- 
Schechter model. 



3 PREVIOUS METHODS 



4 DIFFICULTIES OF TREES 



The development of techniques for constructing Monte Carlo 
r ealizations of the mergin g history of dark matter halos 
(|Kauffmann & White 1993^ |Cole fc Kaiser 1988j |Cole 1991 



Laccy fc Cole 199$ using the extended Press-Schechter for- 



malism has allowed a great deal of progress to be made in 
the use of semi-anal ytic methods for studying ga laxy forma- 
tion and evolution. Kauffmann & White (1993) developed 



a method for constructing merger trees which addresses the 
problem of simultaneously reproducing the average number 
of halos given by Eqn. || and imposing the constraint that 
the mass of a halo be equal to the sum of the masses of its 
progenitors at every stage. To do this they impose a grid in 
mass and redshift. For the first step in redshift, they then 
create a list of halos, where the number of halos with mass 
Mi is given by JV(M<) = iVens AN/AM(Mi) AM,, rounded 
to the nearest integer. iVens is some large number of ensem- 
bles, typically N ens ~ 100. The progenitors are randomly 
assigned to ensembles, starting with the largest and working 
in order of decreasing mass. The probability of assignment 
is proportional to the amount of remaining free mass, with 
the constraint that the total mass of the progenitors cannot 
exceed the mass of the parent^]. This process is repeated for 
all the steps in the redshift grid. 

This algorithm is guaranteed to exactly reproduce the 
mean number of halos of each mass at each redshift for the 
set of ensembles. However, it is possible to encounter a sit- 
uation in which the next halo does not fit into any of the 
ensembles. In order to get around this problem, mass con- 
servation is only enforced in an approximate way (G. Kauff- 
mann, private communication). There are also some practi- 
cal drawbacks to this method. It is necessary to generate a 
large number of ensembles and store them, which is some- 
what inconvenient. An arbitrary grid in halo mass and red- 
shift must be imposed. Also, because the function dN/dM 
is very sharply peaked around Mo for small Afo or small 



redshift intervals Az, the algorithm as described in Kauff- 



mann & White (1993) is sensitive to the binning used and 



is prone to numerical problems. The algorithm breaks down 
for certain choices of power spectrum. Finally, although the 
mean of the distribution is reproduced by construction, the 
partitioning of the halos into individual ensembles is ad-hoc 
and may or may not reproduce the higher moments of the 
distribution. 

A different appro ach, referred to as th e "B lock model" , 
has been proposed by dole & Kaiser (1985) and Cole (1991). 
A major drawback with this approach is that the halo masses 
always grow in discrete steps of factors of two. This is prob- 
lematic for the purpose of semi-analytic galaxy formation 



* Note that in our terminology, the "parent" is younger than its 
"progenitors" , because we always work backwards in time. 



The first choice we must make if we want to construct a 
merger tree is whether to start from high redshift and merge 
together small clumps until the desired redshift is reached 
(as in an N-body simulation or, presumably, in the real uni- 
verse), or, to start from the present day and work backwards 
in time, "disintegrating" the halos into their progenitors like 
a film run backwards. The extended Press-Schechter formal- 
ism provides expressions applicable to both situations. An 
important consideration is that we would like to eventually 
anchor our approach using the observed z = properties of 
galaxies, which are presumably the most secure. In addition, 
because we lack any information about spatial correlations, 
if we go forward in time, we do not know which small clumps 
to combine with which. It may be possible to get around this 
problem somehow, but for now we pursue the "disintegra- 
tion" approach, in which we postulate the existence of a 
"parent" halo of a given mass Mo at redshift zo and break 
it into its "progenitors" working backwards in time. 

It might appear a simple matter to construct a merger 
tree by simply picking the masses of the progenitors of our 
parent halo at some earlier redshift z\ from the expression 
for dP/dM(Mi,Zi\Mo,Zo) given by Eqn. || above, then re- 
peating this process starting from each progenitor in turn 
for the next step back in time. Two difficulties immediately 
arise in implementing this approach. First, from inspection 
of Eqn. ^, the number of halos clearly diverges as the mass 
goes to zero. However, note that the mass contained in small 
halos (Eqn. does not diverge as M — > 0. In order to pick 
masses from the number weighted probability function nu- 
merically, it is necessary to introduce a cutoff mass, or ef- 
fective mass resolution, Mi. 

The second problem is that the progenitor masses must 
simultaneously be drawn from the distribution dP/dM(M) 
and add up to the mass of the parent, Mo. The problem 
is that dP/dM is just the average number of halos that 
one can make out of the mass Mof(M) dAl. What we really 
want is the joint probability function for the set of progen- 
itors {Mi,---,M n }, dP n /dM({Mi,---,M n },Zi | M ,z ) 
with any value of n < Mo /Mi An obvious problem with 
the use of the single halo probability rather than the joint 
probability is that there is no guarantee that we will not 
at some stage pick a progenitor that does not "fit" in the 
halo: i. e. M > Mo — ^ . Mi where Mi are the masses of 
all the previously picked progenitors. In addition, since the 
expression P(M) only gives the probability that there was 
a progenitor of mass M at an earlier time, we do not know 



t Prom now on we will drop the differential notation dP/dM and 
refer to the probability given by Eqn. ^|as simply P(M). We will 
also frequently drop the explicit dependence on the redshift and 
the parent mass Mo where this is unambiguous. 
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a priori how many progenitors were present at the redshift 
z\. 

Having noted these points, we can write down some ba- 
sic requirements for our merger tree construction algorithm: 

(i) The procedure must account for the mass contained 
in halos below the mass resolution Mi, which we will re- 
fer to as the "accreted mass". However the results must be 
independent of the value of Mi. 

(ii) The procedure should treat all progenitors equally, 
independently of the sequence in which they are chosen. 

(iii) The procedure must simultaneously reproduce the 
distribution of the number of progenitors and their masses 
while conserving mass. 

(iv) The algorithm should be numerically robust and 
must be possible to implement in a computationally efficient 
and convenient way. 

We now demonstrate some problems that arise in sev- 
eral seemingly straightforward approaches to building the 
trees. The bold solid line in Figure |l] shows the predic- 
tion of the extended Press-Schechter model for the quan- 
tity AN /AM, the number of progenitors with mass M for a 
parent halo with Mo = 5 Mi after a single step in redshift 
from zq = to z\ = 0.2 (also called the conditional mass 
function). This is the quantity that a successful merging tree 
method must reproduce. In this figure and hereafter unless 
otherwise noted, all masses are given in units of the mass 
resolution Mi. Here we have used Mi = 1.0 X 1O 1O M , but 
the results are independent of this value. 

In our first attempted method, the progenitor masses 
are chosen from Eqn. [B| until the mass reservoir Mo is ex- 
hausted. The probability is set to zero for M < Mi. We have 
tried two ways of addressing the problem of mass "overflow" 
described above. One approach is to choose progenitors until 
the total mass exceeds Mo, and then to truncate the mass 
of the last progenitor. We will refer to this as the "Naive 
Method with Truncation". Another approach is to impose 
an upper mass cutoff, so that the probability is effectively 
set to zero for any values of M that exceed the available 
mass. In effect we then choose the z-th progenitor from the 
modified distribution 

i 

Pi(Mi) = P(M l )0(M l - Mi)Q(M - M i) > ( 6 ) 

3=1 

where P(M) is the original probability function from Eqn. ^ 
renormalized for its new range. No attempt is made to com- 
pensate for the contribution of masses below M; . We refer to 
this as the "Naive Method with Cutoff" . The results for the 
quantity AN /AM (averaged over many ensembles) is shown 
by the histograms marked a (Cutoff) and b (Truncation) in 
Figure |l| 

Both procedures clearly fail to correctly reproduce the 
conditional mass function predicted by the extended Press- 
Schechter model (Eqn. ^) . Although the introduction of the 
upper mass cutoff to prevent choosing masses that are too 
large to fit in the current ensemble is somewhat more elegant 
than the brute force truncation, we see that it produces a 
shift towards smaller masses which leads to large excess of 
small mass halos. It should be noted that this figure shows 
only one step in redshift. The excess multiplies with each 
subsequent step in redshift, so that even a relatively small 
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Figure 1. The number of progenitors with mass m p (condi- 
tional mass function) for a single step in redshift (20 = to 
z\ = 0.2). The bold solid line is the prediction of the extended 
Press-Schechter theory (Eqn. ^). The histograms were obtained 
by picking masses from the distribution P(M,z\ j Mo,zo) until 
the parent mass Mo was exhausted. The dotted histogram (a) 
shows the results using the Naive Method with Cutoff, and the 
short-dashed histogram (b) shows the results of using the Naive 
Method with Truncation (see text). The long-dashed histogram 

(c) uses Accretion Model Method 1, and the dot-dashed histogram 

(d) uses Accretion Model Method 2 (see text). All masses are in 
units of the minimum progenitor mass M; . 

discrepancy quickly becomes very serious. It should also be 
noted that these problems are most pronounced when Mo « 
Mi (as in the case shown in the figure). If Mo Mi, the 
disagreement is not as bad. 



5 THE ACCRETION MODEL 

We must introduce a self-consistent way of treating progeni- 
tors above and below the cutoff mass, Mi . It should be noted 
that although the contribution of masses M < Mi may be 
negligible for halos with Mo Mi , as the mass of the parent 
approaches the resolution limit it necessarily becomes a sig- 
nificant fraction of the total progenitor mass. Because every 
halo, regardless of its size, must be broken down into smaller 
and smaller pieces until all of the pieces fall below the mass 
resolution (this is what makes the process finite) , the correct 
treatment of small halos is crucial for reconstructing the for- 
mation history of halos of all masses out to arbitrarily high 
redshift. 

The basic idea behind this approach is never to treat 
any mass below Mi in terms of progenitor number, but 
rather to find a complimentary description for it as accreted 
diffuse matter. We now introduce an arbitrary distinction in 
terminology to reflect this division. Let progenitors by defi- 
nition have mass greater than the fixed mass resolution Mi . 
The aggregate contribution of all halos with M < Mi will 
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be referred to as accreted mass. A fully rigorous procedure 
should use the joint probability for progenitor number (not 
mass) above Mi, and accreted mass (not number of halos) 
below this mass scale. 

We can now define a few more useful quantities. Given 
the mass of the parent halo Mo and the redshift step 
zo — > z\, the average number of progenitors (recalling our 
definition above), is: 



M 



N={N P (M\M )) = I dM^-P(M, Zl \M ,z ). (7) 

We can also calculate the average fraction of Mo that dwelt 
in the form of progenitor halos of mass M > Mr. 



AM P(M, Zl \M ,z ) 



(8) 



and the complimentary quantity for the average fraction of 
Mo that came from "accreted" mass, / acc = 1 — f p . 

Before we proceed, we would like to warn the reader 
that the contents of the remainder of this section are rather 
detailed and probably only of interest to the specialist. The 
formalism developed in the rest of this section is not used 
directly in the successful method that we will eventually 
derive. We urge the impatient reader to skip directly to 96J. 

From the above predictions we can try to evaluate 
what went wrong with our previous procedures (the "Naive" 
methods). The predicted average number of progenitors for 
the case considered in Figure [l] (Mo = 5 Mi, zo = 0, 
zi = 0.2) is N = 1.14. The actual average number for one 
hundred Monte Carlo realizations using the Naive Method 
with Truncation is N = 2.1, and for the Naive Method with 
Cutoff it is N = 2.3. Thus we see that the mean of the dis- 
tribution is shifted towards larger numbers of small mass 
halos. This motivates our goal in this section, which is to 
find the probability function Vn of having TV progenitors at 
redshift zi given the mass Mo at a later redshift zo, with the 
imposed cutoff of M ; . 

Given Mo, zo, and Z\, let Mi :P be the mass of a progen- 
itor, where by definition Mi iP > Mi- However, during the 
time interval Az = z\ — zo, this progenitor accretes a mass 
Mi :acc due to merging with small halos of mass M < Mi 
which are not counted as progenitors. Therefore its effective 
contribution to Mo is Mi = Ml iP + Mi jacc . We now define 
a modified probability function 



P(Mi|M )= / dMi >P P(Mi, p |Mo)P aC c(Mi, acc |Mi >; 



.(9) 



The weighting function P acc is proportional to the prob- 
ability for Mi >p to accrete a mass Mi jacc during the 
specified redshift interval. This probability is not simply 
P(Mi, acc , Z \\Mo, zo), because M acc will in general be com- 
prised of many small halos. It should be noted that for the 
same reason, M acc is not necessarily smaller than Mi. We 
return to the determination of the function P acc in a mo- 
ment. 

The probability for having one and only one progenitor 
(with additional accreted mass) given a halo of mass Mo is 

Pi=P(M jM ). (10) 

The probability for exactly two progenitors is 



Vi= I dMiP 2 (Mi|M ), 

J Mi 

where 

P 2 (Mi|M ) = Pi(Mi|Mo)Pi(M - Mi | M ) 



(11) 



(12) 



The generalization to TV progenitors Vn is obtained recur- 
sively via 

,Mo-(N-l)M, 

Vn = dM 1 P N (M 1 ,...,M N \M ) (13) 

Jm, 

rMo-iN-^M, 

dMiPi(Mi|M ) 

'M t 

3-1 

r M -J2 Mi-(N-j-l)M, 

*=* dM j P 1 {M j \M ) 

1 

iV-2 

p*0- J2 (M i )-(Af-2)M ! 

i=1 dMjv-i 

IM, 

]V-1 

xPi(Mjv-i \M ) Pi(M - ( M i)\ M °) ■ (I 4 ) 
1=1 

The probability of having no progenitors of mass bigger than 
Ali, Vo, is evaluated at the end by the requirement of 



(15) 



where iV is sufficiently large that Vn — > 0. The average 
number of progenitors must satisfy 



(N p ) = ^2nV„, 



(16) 



and may be compared with the independent prediction of 
Eqn.g 



5.1 The Accretion Probability 

The accretion weighting function P acc (M acc | Al p ;zi,zo) is 
an important missing ingredient in these expressions. It is 
proportional to the probability for a progenitor with mass 
M p to accrete a mass M acc during the redshift interval 
Az — z\ — zo- It should reflect our expectation that it is very 
unlikely for a small mass halo to accrete a very large amount 
of mass, as this would require the simultaneous merging of 
a very large number of halos with M < A/;. Similarly, a 
large halo will be unlikely to accrete a very small amount of 
mass, because its cross section for merging is large. In this 
section we incorporate these qualitative expectations into a 
reasonable guess for the accretion probability function P acc . 

Consider a halo with a mass M at a redshift z^ . F rom 
the spherical collapse model (e.g. White fc Frenk 1991 ), we 
expect the virial mass to increase due to the infall of previ- 
ously uncollapsed material. The mass at a later time corre- 
sponding to a redshift Z2 < z\ is 



(17) 



where the accretion rate from the spherical infall model is: 
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dM,„, v c 3 



(18) 



where V" c is the circular velocity of the halo (the last for- 
mula is strictly true only in a universe with no cosmological 
constant, but is a good approximation even if A 7^ 0). This 
change in mass includes mergers with halos of all masses. 
We still need to estimate how much of the mass change 
AM = Mi — Mi is due to mergers with halos with mass 
less than the resolution limit, i.e. "accretion". To do this we 
use the expression for the mean fraction of accreted mass 
(/acc = 1 — fp, starting from the mass M2 and going back 
in time from 22 to z\. In this way we estimate the average 
mass accreted by M\ to be M acc = / acc (M2) M2. 

For higher moments of the distribution function 
P a cc(M acc |M p ; zi,zo) we shall assume that the accretion is 
mainly due to the infall of blobs (i.e. halos of mass < M;) 
with typical mass Mb <C Mi. For the average accretion of 
Eq. O we expect 



(19) 



For any total accreted mass regardless of the value of Mb, 
the number of blobs is proportional to the average accreted 
mass. If the number of the blobs is Poisson distributed, then 
the second moment of the accreted mass distribution func- 
tion is proportional to the total mass accreted. In the limit 
Nt 3> 1 the distribution should approach a Gaussian dis- 
tribution. We have already argued that the mean of this 
distribution should be A/ acc . We can make a rough guess for 
the width of the distribution, <r acc = /3(AM — M acc ) where 
AM represents the mass change predicted by the spherical 
infall model. This should be an upper limit on the accreted 
mass. Due to the uncertainties involved in the derivation of 
this expression, the parameter (3 is left free and can be tuned 
as needed (we used /3 = 2 for the results presented here) . We 
now have a reasonable guess for the functional form of the 
accretion probability function for a progenitor of mass M p 
over the redshift interval Z\ to zo: 



(M acc - M acc ) 2 

^ u acc 



.(20) 



Pacc(M acc |M p ; zx,Zo) oc — ^ exp 

Although this expression is admittedly ad hoc, one can see 
that it contains the correct qualitative behavior. The mean 
accreted mass increases with the progenitor mass and with 
Az as expected. 



5.2 Merger Trees with the Accretion Model 

We can now imagine a new approach for constructing the 
merger trees which addresses the two main sources of the 
problems in the previous approach: the failure to account 
for accreted mass and the incorrect distribution of the num- 
ber of progenitors. Given the probability function Vn for 
each Mo and time-step, we pick the number of progenitors 
from this distribution. We assign a mass to each of these 
progenitors from the distribution Pi(M) (Eqn. ^) as before. 
The O function prevents us from choosing a mass larger 
than the available mass at any stage. The accreted mass is 
automatically obtained from the residue of this procedure. 
The results of this algorithm (which we will call Accretion 
Model Method 1) are shown by the histogram c in Figure [j] 
We see that the results have improved dramatically from 



the Naive Method with Truncation where the upper mass 
cutoff was also used but the number of progenitors was not 
specified. Also the average number of progenitors is now 
N — 1.08, in much better agreement with the expected 
value of N = 1.14. However there is still an inconsistency 
in this procedure which leads to the remaining discrepancy. 
We have still assigned the mass to the progenitors based 
on the single halo probability function P(M), when, as we 
have argued before, we really should have used the joint 
probability function P n ({Mi, • • • , M n }). This can in princi- 
ple be calculated using the same approach that we used to 
obtain the Vn function. For example the joint probability 
for N = 2, taking into account the accretion weighting, is 
just the integrand of the expression for Vi'- 

P 2 ({Mi,M 2 }) = P(Mi)P(M 2 ) Q(Mo - Mi - M 2 ) 

x / dM acc P acc (M acc I Mi) 

Jo 

xP aC c(M - Mi - M 2 - M acc I M 2 ) . 

This expression may be generalized to the Af-halo joint prob- 
ability as before. 

We now pick the number of progenitors from Vn and 
assign the masses from the joint probability function Eqn. ^l] 
(Accretion Model Method 2). This approach goes a long way 
towards curing the problems we have noted, as we see from 
histogram d in Figure ^. However the shape of the mass func- 
tion is not quite right. We attribute this to the inconsistency 
introduced by our ad-hoc accretion probability weighting. 
We find that changing the form of this function significantly 
affects the results obtained for the mass function. We ob- 
tain better results for a lognormal distribution than for the 
Gaussian distribution used here. If we could somehow ob- 
tain some external information on the form of the accretion 
weighting, for example from N-body simulations, it might be 
possible to produce a successfully working method. However, 
this scheme is also rather cumbersome and computationally 
expensive. The calculation of Vn involves the computation 
of (N — 1)! integrals, and must be repeated for every parent 
halo mass, Mo and redshift interval Az. The joint proba- 
bility function for the i-th progenitor will depend on the 
masses of the previously chosen progenitors and thus must 
be recalculated at each stage. This quickly becomes pro- 
hibitively computationally expensive when large numbers of 
progenitors are allowed. 

We can address the second problem by reducing the 
time-step or redshift interval Az. As Az is decreased, the 
form of Vn steepens and becomes peaked at smaller N. For 
a small enough choice of Az, Vn>s — > 0. We will refer to 
this condition as the two-progenitor limit. The idea is not 
to allow any processes that involve more than a single bi- 
furcation. This demand allows one to calculate only three 
functions for each time step (Vo, Pi, & P2), by using the 
analytic model of Eqn. [l4|. At each stage, the time-step will 
now depend on the parent mass Mo. The larger the halo, 
the smaller the time-step necessary to satisfy this condition. 

Although going to the two-progenitor limit might make 
the procedure computationally feasible, for the moment our 
lack of knowledge about the accretion probability weighting, 
and the sensitivity of the results to this function, lead us 
to abandon this approach. Perhaps the formalism we have 
developed here, and the simple approach we have presented 
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for modeling the joint probability function, can be refined 
in the future. However, we do not pursue it in this paper. 



6 A PRACTICAL SOLUTION 

6.1 Binary Merger Trees Without Accretion 

In the absence of external information about the behavior 
of the accreted mass component, we are forced to treat it 
within our Monte Carlo procedure. As we have discussed, 
we do not want to use the number-weighted probability for 
masses below Mi because of the divergence of this expres- 
sion at small masses. However, the mass-weighted probabil- 
ity, Eqn. |i| does not diverge. If we pick a mass Mi from the 
mass-weighted expression /(Mi, 21 | Mo, 20), this is equiva- 
lent to discovering that a single trajectory, or particle, from 
the parent halo Mo was in a halo with mass Mi at zi. Once 
again, this is the single trajectory probability and if we con- 
tinue to select masses they will not in general fit together in 
any sensible combination that can lead to Mo- We attempt 
to evade this problem by choosing a very small time step 
and so going to the two progenitor limit, as before. It is con- 
venient to use to = S c (z) = S Ct o/D(z) as our time variable 
and S(M) = a 2 (M) as our mass variable as in LC93. These 
can be translated back to redshift and mass by inversion of 
the appropriate expressions. The probability for a step AS 
in a time step Auj is (LC93; Eqn. 2.29) 



P(AS, Au)dAS = 



Auj 



2T(AS)3/: 



■ exp 



IAuj} 



2AS 



dAS (22) 



If we make a change in variables x = Aui/(2\/ AS) this be- 
comes a Gaussian distribution in x with zero mean and unit 
variance. We can see from this expression that if we choose 
the timestep such that 



AM C 



(23) 



where AM C <C Mo, then a step AS corresponding to a 
change in mass AM larger than the mass resolution be- 
comes a 2cr event. We must choose this timestep carefully 
— if it is too big, then the two progenitor approximation 
will break down badly. If it is too small, then the results be- 
come dominated by numerical noise. The above expression 
is approximate, but provides a rule of thumb. Note that 
it scales with Mo through the differential dS/dM(Mo), so 
larger parent halos will require smaller time steps. 

In the simplest version of this algorithm, we start from a 
parent halo with mass A'Iq at 20 and obtain the timestep Auj 
from Eqn. |2^. We work backwards in time from this point. 
We choose a Gaussian random variable with unit variance 
and translate this to a step AS using the transformation 
mentioned above. The new halo mass at the earlier time 
t(ui + Auj) is then M(S + AS). LC93 argue that for a small 
enough timestep, all mergers may be treated as binary. This 
makes the process very simple — at each stage we break 
the halo into two pieces with mass M and AM = Mo — M 
where M is chosen from the probability function f(M). If 
the pr ogenitor obtained in this way is larger than Mi, we 



treat it as the next parent and repeat the procedure. If it is 



smaller than Mi, then we treat it as accreted mass and do 
not follow its history. This is essentially the same algorithm 




Figure 2. The number of progenitors with mass M for a halo 
with initial mass Mq = 500, at various redshifts as shown on the 
figure. The solid lines are the predictions of the extended Press- 
Schechter theory. The histograms show the results for merger 
trees constructed using the Binary Tree Method without accre- 



tion, w hich is essentially the algorithm proposed by Lacey & Cole 



(1993). The trees have an excess of halos compared to the Press- 



Schechter model, and the discrepancy increases with redshift. 



proposed in LC93 at the top of page 641, a nd is similar to 
a gen eraliz ed version o f the block model of Cole & Kaiser 
(1988| ) and |Cole (199l[ ) 



This approach has several advantages. It is simple and 
may be coded recursively in a few lines. Because it mainly 
involves picking Gaussian random deviates, it is also very 
fast. Rather than being imposed on an artificial grid in red- 
shift like previous methods, it reflects the intrinsic merg- 
ing timescales of halos of different mass contained in the 
extended Press- Schechter theory. Unfortunately, the mass 
function of halos obtained in this way begins to develop an 
excess of halos compared to the extended Press-Schechter 
model. This excess becomes more and more severe as the 
number of steps increases. We show this for the conditional 
mass function for a halo with an initial mass Mo = 500A/; 
in Figure ^. This problem becomes quite serious when we 
combine the merger histories of a grid of halos, with the ap- 
propriate Press-Schechter weighting for the parent at z = 0, 
to obtain the total mass function, shown in Figure [| The 
mass function reconstructed in this way should yield exact 
agreement with the original Press-Schechter expression for 
the universal mass function. With this method, the number 
of halos is overpredicted by almost an order of magnitude 
by a redshift of 2 ~ 6. This is especially troublesome as the 
Press-Schechter model already gives an ~30- 50% excess of 



smal l halos compared t o N-body simulations ( Lacey fc Cole 
1994; Gross et al. 1998). The pro blem that we demons trate 



here may explain why LC93 and Lacey & Cole (1994) find 
that a similar Monte Carlo method leads to halo formation 
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10 11 12 13 



log M (M e ) 

Figure 3. The total mass function obtained from the Binary 
Tree Method without accretion. All mergers are assumed to in- 
volve exactly two halos. The total mass function is obtained by 
combining a grid of halos from lMj to5xl0 4 Af; (l.OxlO lo M to 
5.0 x 10 15 M Q ) weighted with the Press-Schechter number density 
at z = 0. The bold line is the prediction of the Press-Schechter 
theory, at z = 0.16, 2.5, 3.6, 5 and 7 from top to bottom. The 
merger trees (histograms) overpredict the number of halos by 
more than an order of magnitude after many steps in redshift. 

times that are 40% higher than the analytic predictions or 
the N-body results. 

We believe that this problem is due to the simplifying 
assumption of binary mergers. Because the merger rate of 
very small halos becomes effectively infinite for CDM-like 
power spectra, non-binary mergers, at least involving small 
halos, cannot be neglected. This statement is complemen- 
tary to our original premise regarding the importance of 
what we have called accreted mass. One might think that 
we have simply not chosen a small enough time step, but 
if this were the case the results should improve steadily as 
we decrease the time step. We do not observe this behavior 
even for extreme reductions in the time step. 

Another way of stating the problem is that we have 
actually violated item 2 in our list in §^j. The first mass 
is chosen from the distribution f(M), but the mass of the 
second progenitor is not, it is just assumed to be whatever 
mass is left over. This means that for every progenitor with 
mass M we always get a progenitor with mass Mo — M. It 
is easy to see that this will lead to inconsistencies with the 
mean distribution function P(M). 

6.2 Binary Trees with Accreted Mass 

We now attempt to cure the problem noted above by re- 
laxing the simplifying assumption of binary mergers from 
the previous subsection. Namely, we postulate that mergers 
can involve at most two progenitors (halos with M > Mi), 
but an arbitrary number of halos with mass less than Mi. 



This of course amounts to allowing for accreted mass. At 
any branching we may have only accreted mass (zero pro- 
genitors), or alternatively one or two progenitors plus ac- 
creted mass. In addition, the progenitor masses must always 
be picked from the probability distribution ,f(M). Leftover 
mass can contribute to accretion but cannot be used for 
progenitors. 

The new recipe is as follows. The algorithm is shown 
in flow-chart form in Figure ^. Given the parent mass Mo 
we compute the time step Auu as before. Using this timestep 
throughout the following steps, we: 

(i) Pick a mass M from the mass-weighted probability 
distribution Eqn. This mass can be anywhere in the 
range < M < Mo- If M < Mi, we count it as accreted 
mass. If M > Mi, we count it as a progenitor. 

(ii) Compute the unallocated mass AM = Mo — M. 

(iii) If the unallocated mass AM is larger than Mi, then 
it may or may not contain a progenitor. To determine this, 
pick another mass M from the distribution, but with the 
restriction M < AM. Depending on its mass, count it as 
accreted mass or a progenitor as before. In either case, sub- 
tract M from the mass reservoir. 

(iv) Repeat this process until either 

• The mass reservoir AM falls below the minimum 
halo mass M; , in which case it must abandon any aspira- 
tions of harboring a real progenitor and must be accreted 
mass, 

• OR we have found a total of two progenitors (M > 
Mi), in which case the remaining mass is considered to 
be accreted mass in accord with our ansatz. 

(v) Each progenitor now becomes a parent, we calculate 
a new time step, and repeat the whole process. 

In the flow-chart, branches leading to the outcome of zero, 
one and two progenitors are labeled Vo, V\, and V2, in con- 
nection with the formalism developed in the previous sec- 
tion. 

Note that this procedure does not strictly fulfill require- 
ment 2 of equal treatment of progenitors regardless of the 
order in which they are picked. This inequality is necessary 
due to the mass conservation requirement. 

The results for the conditional mass function of a single 
halo are shown in Figure The discrepancy is now in the 
opposite direction: the number of halos is underpredicted 
relative to the extended Press-Schechter prediction. Appar- 
ently this procedure now overestimates the accreted mass. 
This is not too surprising since we allowed large amounts of 
mass to be designated as accreted mass simply because two 
progenitors had already been found. It appears that even 
in the limit of small time steps and for large mass halos, 
mergers between more than two halos cannot be neglected. 
In some ways this is not surprising either, because after all 
the division into M > Mi and M < Mi is arbitrary and has 
no physical basis. 

6.3 The Successful Method: TV-Branch Trees With 
Accretion 

It is trivial to generalize our previous recipe to allow an un- 
restricted number of progenitors. We now continue picking 
progenitor masses until the unallocated mass AM is less 
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Figure 4. A flow-chart for one redshift step of the disintegration of a halo in a merging tree. The end points (V\,'P2, etc.) lead to 
identical flo w-ch arts for the subsequent timestep, with a new parent mass Mo. A detailed discussion of the algorithm and this flow-chart 
is given in 36.2. 



than Mi. This is indicated on the flow-chart by the dashed 
line labeled "To 7V' . Note that the total accreted mass can 
still exceed Mi because some of the attempts to pick pro- 
genitors yield halos with M < Mi and contribute to the 
accreted mass. We still pick the time-step so that the num- 
ber of progenitors cannot get too large (we find that we 
never exceed ten progenitors per time step even for cluster 
mass halos (Mo = 5 x 10 4 Mi)). We find a good compromise 
between efficiency and accuracy if we introduce an addi- 
tional scaling in the expression for Alj from Eqn. of the 
form 6 + alog 10 (M/M;), where we have used the parameters 
a — 0.3 and b — 0.8 for the results shown here. This optimal 
scaling would change for a different power spectrum shape. 

This recipe gives good results for the conditional mass 
function for parent halos with a wide range of masses. 
We show this in Figure ^| for parent halos with Mo = 
5 Mi, 500 Mi and 5 x 10 4 M;. The agreement is poorest for 
halos with Mo < 10 M; . If we require strict mass conser- 



vation, this is an unavoidable problem due to the shape of 
the conditional mass function for halos of this size. This 
should be kept in mind when setting the value of Mi — it 
should be chosen such that only objects larger than ~ 10M; 
correspond to observable galaxies. We also check the mass 
weighted quantity f p , the fraction of mass in progenitors as 
a function of redshift. This is shown in Figure [jj. This quan- 
tity shows good agreement for the smallest halo, Mo = 5 Mi, 
which shows that the accreted mass is being treated prop- 
erly, so that we do not need to worry about the less than 
perfect agreement in the conditional mass function, as long 
as the condition on Mi mentioned above is satisfied. Note 
that the agreement of the mass function can be improved 
by adjusting the time-step, but at the expense of f p . We 
adjust the time-step to achieve the best possible agreement 
for both the number and mass weighted quantities, over the 
entire mass range. We compute the universal mass function 
from the weighted grid of merging histories constructed us- 
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Figure 5. The number of progenitors with mass M, for a halo 
with Mo = 500, using the Binary Tree Method with accretion. 
Here mergers may involve only two halos with mass greater than 
M[, but an indefinite number of mergers with halos with mass 
smaller than M; (accretion). As usual the solid line is the ex- 
tended Press-Schechter prediction. The trees (histograms) now 
underproduce halos compared to the Press-Schechter model. 



ing our new scheme, and plot this in Figure g. We now find 
very good agreement with the prediction of the standard 
Press-Schechter theory. 

We therefore conclude that although this method is not 
rigorous, it produces acceptable agreement with the mean 
quantities that we can check with the Press-Schechter model. 



7 THE NUMBER-MASS DISTRIBUTION OF 
PROGENITORS 

We have now developed a convenient and efficient method 
for constructing merger trees. The averages derived from 
an ensemble of these trees agree with the important mean 
quantities predicted by the extended Press-Schechter the- 
ory. However, part of the motivation for developing this 
new method was to ensure that the ensemble obeys the true 
joint probability distribution. We have not yet shown this 
to be the case, and we have already mentioned the lack of 
any information from the standard Press-Schechter formal- 
ism which would allow us to evaluate the Monte Carlo trees 
against analytic predictions. 

We do have the predictions of the model we developed in 
section ^, but we do not trust them for reasons discussed in 
that section. However, out of curiosity we compare the pre- 
dictions of the semi-analytic accretion model for the proba- 
bility distribution of the number of progenitors, Vn, with the 
results of the Monte Carlo merger trees. The integrals were 
performed numerically using a recursive, adaptive stepsize 
Runga-Kutta algorithm. Figure ^ shows the Vn distribu- 
tions for halos with mass Mo = 5 and Mq = 50 at redshifts 
of 0.1, 0.2, 1, and 3. The results agree rather well for redshift 




1 10 100 1000 10" 1 10 100 1000 10" 



Figure 6. The number of progenitors with mass M, for a halo 
with initial mass Mo, using the TV-Branch Tree Method with ac- 
cretion. This is the same as the binary tree method with accretion 
except that an arbitrary number of progenitors is allowed at each 
branching. The solid line is the extended Press-Schechter predic- 
tion, (a) M = 5 Mi (b) M = 500 Mi (c) M = 5 X 10 4 Mi. The 
merger trees (histograms) are in reasonably good agreement with 
the extended Press-Schechter model. 
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Figure 7. The average fraction of the original mass Mo contained 
in progenitors (halos with M > Mi) at a redshift z. The solid line 
shows the prediction of the Press-Schechter model. The square 
symbols show the average given by the merger trees (TV-Branch 
Tree Method with accretion) Error bars show the standard devi- 
ation (lcr) over many ensembles. 



Figure 9. The probability distribution of the number of progen- 
itors, for a parent halo with mass Mo at several redshifts. The 
histogram shows the results from the Monte Carlo merging trees 
(N-Branch Tree Method with accretion). The stars show the re- 
sults of the semi-analytic model developed in section ^ Left panel: 
M = 5 Mi . Right Panel: M = 50 M t 




10 11 12 13 14 

log M (MJ 



Figure 8. The total mass function obtained from the TV-Branch 
Tree Method with accretion (the successful method). The solid 
lines show the prediction of the Press-Schechter theory, at z = 
0.16, 2.5, 3.6, 5 and 7 from top to bottom. Broken histograms 
show the mass function from the merger trees, constructed as 
described in Fig. H. 



Figure 10. The evolution with redshift of the two-dimensional 
distribution of the number of progenitors and their masses for a 
halo of mass Mo, obtained from Monte Carlo realizations of the 
TV-Branch Merger Trees with accretion, (a) Mo = 5 M; (b) Mo = 
50 Mi (c) Mo = 500 Mi. These figures are available by anonymous 
ftp from ftp.ucolick.org (pub/outgoing/rachel). 

steps of Az = 0.2 or smaller. For larger steps in redshift, the 
semi-analytic results do not agree as well, probably due to 
the breakdown of the accretion model. 

Figure [ji] demonstrates the importance of taking into 
account the joint distribution in the progenitor number - 
progenitor mass space. The figure shows strong correlations 
between the two variables. Some of the correlations are ob- 
vious: in Figure [ic|a, we show the distribution for a par- 
ent halo of Mo = 5 Mi. In the first redshift step (20 = 0, 
Zi — 0.2) the highest probability is obtained for having a 
single progenitor, and this progenitor naturally contains a 
large fraction of the mass Mo. As we progress in redshift, 
this correspondence is not preserved. Accreted mass starts 
to be more and more significant, and the unseen accreted 
mass complements low mass progenitors so that the sum 
may reach Mo. At the formation epoch of Mo, all Vn are 
populated, and in earlier stages most of the halos go below 
Mi, when the dominant process is of single progenitors that 
accumulate accreted mass. It is interesting to notice that 
the "formation epoch" (the earliest time when the largest 
progenitor has mass greater than Mo/2) is not dominated 
by mergers of equal mass halos, but rather approached via 
slow accretion. Thisgeneral picture also remains valid for 
higher Mo (Figure [lOp and c): an increase in the number 
of progenitors occurs towards an intermediate redshift, and 
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it then declines towards containing most of the mass in the 
accreted component. However it should be noted that the 
highest mass considered here would be comparable to the 
halo of a Milky Way sized galaxy (5 x 1O 12 M ). For much 
larger mass halos (comparable to group or cluster mass ha- 
los), accretion is less important relative to the aggregation 
of roughly equal mass progenitors. 

The probability for the number of progenitors spans 
substantial parts of its permitted range even for Mo = 500. 
It is therefore clear why an infinitesimal step is needed for 
the two-progenitor scheme to work. As soon as we consider 
a finite timestep, the probability for Vn>i is no longer neg- 
ligible. 

The details of the redshift sequence represent the char- 
acteristics of the specific power spectrum and cosmology we 
used for the Monte-Carlo realization. The formation time as 
a function of halo mass and redshift are determined by the 
cosmology and the power spectrum. The qualitative trend 
of this sequence, however, is similar for all cosmologies and 
all hierarchical power spectra. 

Figure suggests that the progenitor mass and number 
distribution functions are an interesting avenue to pursue in 
the study of structure formation via merger trees. More im- 
portantly, it points at the existing interplay between the 
accreted mass, the progenitor mass and the number of pro- 
genitors; an interplay in which none of the three can be 
treated separately from the other. 



aries of the existing theory, and so must be examined by 
comparisons with N-body simulations. However, the simu- 
lations have their own problems and complications, such as 
the limitations of mass and spatial resolution and the ambi- 
guities of defining halos, so they should not be regarded as 
necessarily representing the absolute truth. In addition the 
agreement between the simulations and the Press-Schechter 
model is only approximate, even for the mean quantities 
such as the mass function. It would therefore be desirable 
to have a reliable theoretical means of addressing this prob- 
lem. We have attempted to reformulate the extended Press- 
Schechter theory to obtain the full probability distribution 
for the number of progenitor halos Vn- Although this model 
gives qualitatively reasonable results for certain cases, some 
ingredients remain ad-hoc. 

For the moment, this leaves us with no recourse but 
t o appeal to N-body simulations. In a companion paper 
flSomerville et al. 1998[ ), we will compare the results we have 
obtained here with numerical simulations. This comparison 
has two goals: (a) to determine the quality of agreement 
of the analytic and Monte Carlo results with the N-body 
simulation results, (b) to study the full distribution of the 
various quantities, and determine whether the Monte Carlo 
method developed here reproduces these results. The merger 
trees will then be used as the framework for the development 
of full semi-analytic galaxy formation mo dels, and used to 
compare with a variety of observations ( [Somerville 199 ' 
Somerville & Primack 1998). 



8 SUMMARY AND CONCLUSIONS 

We have presented a new method for constructing the merg- 
ing history of dark matter halos in a semi-analytic way. We 
have highlighted the need to impose an arbitrary mass cutoff 
for practical reasons, which leads us to distinguish between 
halos above and below this threshold as "progenitors" or 
"accreted mass" , respectively. The scheme we have proposed 
and implemented treats accreted mass and progenitors in 
a self consistent way, and produces good agreement with 
the average quantities predicted by the underlying Press- 
Schechter theory, such as the conditional and universal mass 
function of halos and the mean mass in progenitors as a func- 
tion of redshift. 

Our method is an im provement on the method proposed 
by Lacey & Cole (1993), which, after many steps in red- 
shift, substantially overproduces halos relative to the Press- 
Schechter mass function. Our work suggests that it is not 
possible to simultaneously conserve mass exactly and retain 
the exact agreement with the conditional mass functio n from 
the extended Press-Schechter model. The method of Kauff- 



mann & White (1993) reproduces the conditional mass func- 



tion exactly and conserves mass approximately. Our method 
conserves mass exactly and reproduces the conditional mass 
function approximately. This seems to be a necessary trade- 
off. Our method does have certain practical advantages: it 
does not require the creation and storage of a large number 
of ensembles, it is numerically robust, it does not require the 
imposition of a grid is mass or redshift, and it will work for 
any power spectrum. 

We have pointed out the necessity of investigating the 
full probability distribution of the number of progenitors 
and their masses. This cannot be tested within the bound- 



ACKNOWLEDGEMENTS 

We would like to thank Gerard Lemson for useful discussions 
and healthy skepticism. We also thank Guinevere Kauff- 
mann and Simon White for their thorough reading of the 
manuscript, an important demonstration of some drawbacks 
in the proposed procedure, and comprehensive discussions. 
We thank Sandy Faber and Joel Primack for detailed com- 
ments on the manuscript and useful insights. RSS acknowl- 
edges support from an NSF GAANN fellowship. This work 
was supported in part by NASA ATP grant NAG5-3061. 



REFERENCES 

Bardeen J. M., Bond J. R., Kaiser N., Szalay A. S., 1986, ApJ, 
304, 15 (BBKS) 

Bond J. R., Cole S., Efstathiou G., Kaiser N., 1991, ApJ, 379, 
440 

Bower R., 1991, MNRAS, 248, 332 
Cole S., 1991, ApJ, 367, 45 

Cole S., Aragon-Salamanca A., Frenk C, Navarro J., Zepf S., 

1994, MNRAS, 271, 781 
Cole S., Kaiser N., 1988, MNRAS, 233, 637 
Efstathiou C, Bond J., White S., 1992, MNRAS, 258, 1 
Efstathiou C, Frenk C. S., White S. D. M., Davis M., 1988, MN- 
RAS, 235, 715 
Eisenstein D., Loeb A., 1996, ApJ, 459, 432 
Gelb J. M., Bertschinger E., 1994, ApJ, 436, 467 
Gross M. A. K., Somerville R. S., Primack J. R., Holtzman J., 

Klypin A. A., 1998, MNRAS, submitted, astro-ph/9712142 
Kauffmann G., 1996, MNRAS, 281, 487 
Kauffmann G., White S., 1993, MNRAS, 261, 921 
Kauffmann G., White S., Guiderdoni B., 1993, MNRAS, 264, 201 



© 0000 RAS, MNRAS 000, 000-000 



How to Plant a Merger Tree 

Lacey C, Cole S., 1993, MNRAS, 262, 627 (LC93) 
Lacey C, Cole S., 1994, MNRAS, 271, 676 
Ma C, 1996, ApJ, 471, 13 

Peebles P. J. E., 1980, The Large-Scale Structure of the Universe. 

Princeton Univ. Press, Princeton, NJ 
Press W., Schechter P., 1974, ApJ, 187, 425 

Somcrville R., 1997, Ph.D. thesis, University of California, Santa 
Cruz, http://www.fiz.huji.ac.il/ rachels/thesis.html 

Somerville R. , Lemson G., Kolatt T., Dckcl A., 1998, in prepara- 
tion 

Somcrville R. S., Primack J. R., 1998, MNRAS, submitted 
White S., Frenk C, 1991, ApJ, 379, 52 



© 0000 RAS, MNRAS 000, 000-000 



